function [V_c,V,Cell] = volume
%计算体积与体心坐标
%   此处显示详细说明
global own nei num d f_c c_s
V_c=zeros(num,3);
V=zeros(num,1);
V_c_=zeros(num,3);
V_num=zeros(num,3);
Cell=cell(num,1);
for i=1:own(1)
    V_c_(own(i+1),:)=V_c_(own(i+1),:)+f_c(i,:);
    V_num(own(i+1),:)=V_num(own(i+1),:)+1;
    Cell(own(i+1))={[cell2mat(Cell(own(i+1))),i]};
end
for i=1:nei(1)
    V_c_(nei(i+1),:)=V_c_(nei(i+1),:)+f_c(i,:);
    V_num(nei(i+1),:)=V_num(nei(i+1),:)+1;
    Cell(nei(i+1))={[cell2mat(Cell(nei(i+1))),i]};
end
V_c_=V_c_./V_num;
for i=1:num
    lface=cell2mat(Cell(1));
    n=size(lface,2);
    V_=zeros(n,1);
    D=zeros(n,3);
    for j=1:n
        D(j,:)=0.75*f_c(lface(j),:)+0.25*V_c_(i,:);
        V_(j)=abs(c_s(lface(j),:)*(f_c(lface(j),:)-V_c_(i,:))'/3);
    end
    V_c(i,:)=V_'*D/sum(V_);
    V(i)=sum(V_);
end
d=zeros(num,num);
for i=1:nei(1)
    d(own(i+1),nei(i+1))=sqrt(sum((V_c(nei(i+1),:)-V_c(own(i+1),:)).^2));
end
d=d+d';
end
